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Abstract 

Let us consider some Coulomb systems of several infinitely massive centers of charge Z and one- 
two electrons: {Z,e), (2Z, e), {3Z,e), (4Z, e), (2Z, e, e), (3Z, e, e). It is shown that the physical, 
integer charges Z = 1,2, ... do not play a distinguished role in total energy and equilibrium distance 
of a system giving no indication to a charge quantization. 

By definition, a critical charge Zcr for a given Coulomb system {nZ, e) or {nZ, e, e) is a charge 
which separates a domain of the existence of bound states from a domain of unbound ones (con- 
tinuum). For all above-mentioned systems critical charges Zc as well as equilibrium geometrical 
configurations are found. For all studied systems there was obtained an indication to a square-root 
singularity at Z = Zcr- 
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I. INTRODUCTION 



One of the basic observations of fundamental physics is the quantization of electric charges 
of elementary particles and nuclei. The electric charges of electron and proton have opposite 
signs and their values coincide, neutron has zero electric charge, any nuclear electric charge 
is equal to proton charge multiplied by integer number. This observation is supported 
experimentally, gets its justification in elementary particle theory and nuclear physics. A 
natural question to ask is there any indication to such a quantization in atomic-molecular 
physics. In classical electrostatics the stable configurations of point charges are absent 
(the Earnshaw's theorem), zero charge looks like as a singular point where the nature of 
interaction changes from repulsion to attraction. Usually, at a singular charge the whole or 
a part of the potential vanishes. In non-relativistic quantum electrodynamics these singular 
charges undoubtedly continue to exist, however, a new phenomenon occurs - there are some 
critical charges which separate the domain of the existence of the bound states from the 
domain of non-cxistencc, although the nature of potential remains unchanged. In some 
cases a system gets bound at a critical charge with polynomially-decaying eigenfunctions 
at large distances unlike standard exponentially-decaying eigenfunctions. To the best of 
our knowledge this phenomenon was observed for the first time for the Helium-like system 
(Zee). It was named as "the level hits (kicks) continuum", or as "the zero-energy state", or 
as "the level on the threshold of continuum". Probably, two simplest examples where such 
a phenomenon occurs arc the Poschl- Teller potential and the Yukawa potential. 

In this paper wc consider a Coulomb system of a number of infinitely massive centers of 
the same charge Z and one- two electrons assuming that the charge Z is a real parameter. 
The main goal of this paper is to explore a question: are integer (physical) charges Z special 
in some sense when the total energy is studied? Another goal is to find the domain(s) in Z 
where the system has at least one bound state. We intend to find the critical charges Zcr 
which separate the domains of existence/non-existence of bound state. The study is made 
in framework of the non-relativistic quantum mechanics. 
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II. GENERALITIES 



Let us consider the Coulomb molecular system which consists of n fixed charges Z and k 
electrons, {nZ, ke). The Hamiltonian which describes this system is written as follows 

-I ^ fy2 n k ^ k ^ 

« = 4E^. + + (1) 

a=l i<j ■> 1=1 a=l a<b 

in a.u., where Rij is the distance between charge centers i and j, is the distance from 
ath electron to «th charge center, Tab is the distance between electrons a and 6, and k = 1,2. 
If Z = 1, the Hamiltonian describes the system of n protons and one-two electrons in 
the Born-Oppenheimer approximation of the zero order (the protons are considered to be 
infinitely-massive). There are three important particular cases. 
(1) Atomic-type case, n = 1. The Hamiltonian ([T]) gets a form 

a=l a=l a<b 

where is the distance between zth electron and the center. At A; = 1 (one-electron case) 
we get a hydrogen-like ion, its spectrum is known 

En{Z) = , (3) 

where is the principal quantum number, = 1, 2, . . .. Discrete spectra is infinite for any 
Z > 0. Critical point is at Zcr = 0. Nature of this critical point is of quite obvious - it is 
a singular point of the differential equation, at Z = the potential vanishes. It is worth 
noting that En{Z) has no singularities at finite Z having the pole of the second order at 
Z = oo. For n > 1 it is evident from physical point of view that for small Z the system 
is unbound but gets bound for sufficiently large Z. Hence, there exists some Z = Zcr- 
It seems established that Zcr{n = 2) ^ 0.91 (see e.g. l|, 121] and references therein) and 
Zcrin = 3) < 2.1 (see e.g. |3!])- Making a rescaling of the Hamiltonian (|2]), r -I we get 
the Hamiltonian in the form 

-M k k -M k ^ 

a=l a=l a<b 

and arrive immediately at the conclusion that the energy of the bound state has the second 
order pole at Z = oo. In general, in the domain [Zcr, oo) the ground state energy E{Z) is 
smooth monotonous function of Z without any indication to a charge quantization. 



(2) One-electron, molecular-type case, k = 1. The Hamiltonian has a form ([T]) without 
the last sum. It is evident that for small Z the system is bound and one critical charge 
coincides with the singular point of the Hamiltonian Z^r = 0, where the potential vanishes. 
For large Z the Coulomb repulsion of charged centers gets larger than the attraction of 
the electron to them and a system definitely becomes unbound. Hence, the second critical 
charge at finite Zcr > must exist. Our goal is to find this critical charge for n = 2,3,4. 

(3) Two-electron, molecular-type case, k = 2. The Hamiltonian ([1]) gets a form 

From physical point of view it seems evident that a system is not bound for large Z as well 
as for Z < 0. Thus, there must exist two critical points: one has be near zero, Z ~ and 
another one has to be at finite Z. None of them is of a type of singularity of the operator 
(|5]). Our goal is to find this critical charge for n = 2, 3. 

It is necessary to introduce a formal definition of the critical charge Z = Z^- for molecular 
system. It is natural to do it in the Born-Oppenheimer approximation of the zero order when 
Z charges are assumed to be fixed. In the case of the existence of a bound state the potential 
curve Etotai = Et{R) has a minimum at finite internuclear distance R = Reg. If the bound 
state is stable the potential energy at infinite intercenter distance is larger at Rgg. Otherwise, 
the bound state can be metastable globally: the system can decay and the potential energy 
at infinite intercenter distance is smaller than at minimum. It implies the existence of the 
maximum on the potential curve at some finite R > Req. In the case of non-existence of 
a bound state the potential curve has no minimum at finite distances. Hence, at critical 
charge Z = Z^r the potential curve has a saddle point at a finite distance R. 
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III. ONE-ELECTRON MOLECULAR SYSTEMS 



A. Two center case (2Z, e) 




FIG. 1: Geometrical setting for (2Z, e) system 



It is well known that at Z = 1 there exists the molecular hydrogenic ion H^, while at 
Z = 2 the molecular helium ion Hcg^ does not exist. Hence, the critical charge Zcr has to 
be in the range 1 < Z < 2. 

In order to calculate the total energy of {2Z, e)-system vs i? as a function of charge 
Z we use the variational method. As a trial function is taken linear superposition of the 



Heitler-London, Hund-Mulliken and Guillemin-Zener functions (see {4]). 
I. The Heitler-London function. 



^1 



-aiZ{ri+r2) 



where ai is variational parameter. It is worth mentioning the potential, for which the 
function \E^i is exact ground state wavefunction, 

reproduces both Coulomb singularities and at ai = 1 even their residues. The parameter 
«! 7^ 1 makes sense ( ant i) screening of the nuclear charges. It is well-known for Z = 1 that 
the Heitler-London function describes small internuclear distances and can give a significant 
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contribution near equilibrium, at i? ~ Req. It mimics a coherent interaction of the electron 
with charged centers. It seems evident it holds for Z ^ 1. 
II. The Hund-MuUiken function. 



^2 = (e""2^^i + e-"2^^2) 

where a2 is variational parameter. It describes incoherent interaction of the electron with 
charged centers. This function gives a significant contribution for large internuclear dis- 
tances. 

In order to interpolate between domains R ~ R^g and R ^ Reg, we use two interpolating 
functions. 

III-l. The Guillemin-Zener function 

It is the simplest non-linear interpolation between -01 and or, saying differently, be- 
tween small and large internuclear distances, 

_ ^^-a3Zri-a4Zr2 _|_ ^-asZr2-a4Zri^ 

where as, 04 are variational parameters. If 

• CKs = 0^4 then -03^ — )■ i/ji 

• 04 = then ips-^ — > 'ip2 
III-2. Linear Interpolation 

IV. Superposition of the two kinds of interpolation 

-04 = As^tjjs^ + As^ijjs^ . (6) 



With such a six-parametric trial function ([6]) IJ] the expected relative accuracy in total 
energy is ~ 10^^, which is confirmed by an independent calculation based on use of highly 
accurate uniform approximation of the ground state eigenfunction [g] (see a discussion be- 
low). The total energy E{Z,R = R^g) is presented at Fig. [2] and the equilibrium distance 
is at Fig. [3l Both curves are smooth without any indication to a special feature at the 
physical charge Z = 1. At some charge (see below) the energy curves for {2Z, e) and Z-atom 
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{Z, e) intersect. This crossing separates the domain of stabihty from metastabihty of the 
system (2Z, e). The equihbrium distance vs. Z is a smooth curve which has a minimum 
Reg = 1.952 a.u. at Z = 0.7924, expectedly, with a decrease of Z it grows to infinity. At 
critical charge Z = 1.439 the equihbrium distance, where the potential curve E = E{R) has 
the saddle point, is equal to 2.985 a.u. 

It is interesting to study the approach of the total energy to the critical charge Z — > Zcr 
from below. In order to do it we use the Puiseux expansion 

oo 

ET{Z) = Y,Cin{Zcr-Zf" , (7) 

n=0 

with the condition that fe„ < Our goal is to find parameters a„ and 6„ of this expansion. 
Restricting the expansion ([7]) to a finite number of terms we make fit of the total energy 
calculated numerically, see Table I. The fit based on data from the domain Z G [1.30, 1.43] 
(20 points) is: 

Et{Z) = -1.8072 + 1.5538 {Z^r - Z) - 0.5719 - Zf^ (8) 
+ 0.1129(Ze^ - Zf + 0.7777(Ze^ - Zf'^ - 0.4086(Zc^ - Zf^ + ... 

where the critical point is 

Z'^^ = 1.439 . (9) 

This behavior indicates that critical point might be a square-root branch point. 

There exists a charge for which a type of the binding of the system (2Z, e) is changed 
from metastable, (2Z, e) — > {Z, e) + Z, to stable, Z = Zcross = 1-237 at Req,cross = 2.184 a.u. 
It corresponds to the crossing of two potential curves on Fig. O If Z < Zcross the system 
(2Z, e) is stable, ii Z > Zcross the system (2Z, e) gets unstable, (2Z, e) — > (Z, e) + Z. Value 
of Zcross "we calculated coincides with one found by Rebane [3]. 

As next we study the behavior of the total energy near the point of crossing, Zcross- From 
the left , Z < Z cross-, we find as the result of the fit that the Puiseux expansion becomes the 
Taylor expansion 

Et = -1.5292 + 1.341 (1.2366 - Z) + 0.08 (1.2366 - Z)^ + . . . 

as well as from the right, Z > Zcross, our data are also fit by the Taylor expansion 

Et = -1.5292 + 1.340 {Z - 1.2366) + 0.05 {Z - 1.2366)^ + . . . 
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z 




Fit 


0.10 


-0.031019 


-0.03071 


0.15 


-0.064596 


-0.06455 


0.20 


-0.107149 


-0.10735 


0.25 


-0.157038 


-0.15725 


0.30 


-0.212917 


-0.21287 


0.35 


-0.273656 


-0.27336 


0.40 


-0.338292 


-0.33838 


1.30 


-1.614220 


-1.61422 


1.32 


-1.641112 


-1.64111 


1.34 


-1.668126 


-1.66813 


1.36 


-1.695327 


-1.69533 


1.38 


-1.722801 


-1.72280 


1.40 


-1.750671 


-1.75067 


1.41 


-1.764813 


-1.76481 


1.42 


-1.779144 


-1.77914 


1.43 


-1.793737 


-1.79373 



TABLE I: Total energy of {2Z, e) in Ry at equilibrium vs Z obtained using 
method ^ compared to the result of the fit ([8]) . 



and in the 



Inside of the accuracy of data used these expansions do coincide. Therefore, we do not see 
an indication to a branch point singularity contrary to the statement in |3|. It is worth 
mentioning that the dependence of Req on Z near Zcross is also very smooth, see Fig. [3l 

Another question to rise is a behavior of the total energy near point Z = which is 
the singular point of the Schrodinger equation. Based on the fit of data from the domain 
Z G [0.1,0.5] (six points, see Table I) we find that the Puiseux expansion ([7]) becomes the 
Taylor expansion 



-3.5258Z^ + 4.8922Z^ - 3.4121Z* + . . . . 



Such a behavior does not provide an indication to a singular nature of the point Z = 0. 
However, the total energy can not be analytically continued to ReZ < 0. 
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FIG. 2: Total energy in Ry of systems (Z, e) (dotted line), (2Z, e) at i? = R^.^ (dashed line), (3Z, e) 
at i? = Req (long-dashed line) and (4Z, e) at = i?eg (solid line) as functions of the charge Z. 
(Z, e) and (2Z,e) curves cross at Z = Zcross = 1.237. Dashed curve ends at Z 
Long-dashed curve ends at Z = Zc? = 0.9539. Dotted curve ends at Z = zfr = 0.736 



Z^cr = 1.439. 



B. (3Z,e) 



Let us consider the electron in the electric field of three static charges Z: (3Z, e). In 
general, these charges form triangle, see Fig. H] as an illustration. Such a system does not 
exist at Z = 1 [8]. Thus, there might exist a critical charge Z < 1 for which the system 
gets bound, it separates the domain of the non-existence from existence of the bound state. 
Evidently, one of such critical charges is at Z = 0, which is the singular point of the 
Schrodinger equation. Another one is at some Z = Z^r < 1 (see [s]). Calculations (see 
below) show that Zcr = 0.9539 with Reg = 4.754 a.u. Thus, the system (3Z, e) exists for 
charges < Z < Zcr always in a form of equilateral triangle, which is the optimal geometrical 
configuration. It was checked that this configuration is always stable with respect to small 
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1.6 



0.4 



z 



0.8 



1.2 



FIG. 3: Equilibrium distance in a.u. of systems (2Z, e) (dashed line), (3Z, e) (long-dashed line) 



and (4Z, e) (solid line) as functions of the charge Z . All curves cross at Z = 0.2977 with 
2.366 a.u. Dashed curve ends at Z = Zcl^ = 1.439. Long-dashed curve ends at Z = Zcr^ 
Dotted curve ends at Z = zif^ = 0.736. 



cross 
■eq 



0.9539. 



deviation. 

In order to calculate the total ener gy E {Z, R) the variational method is used. We employ 



the physics-inspired trial functions 



ll| taking afterwards their linear superposition, 



trial 



(10) 



here Aj are linear parameters. Each function ij)^^^ is chosen in such a way to describe a 
certain physical situation of the system. In general, il)^^"^ has the form of a symmetrized 
product of three Is-Coulomb orbitals (Slater functions) 



_ g-airi-a2»'2— oaf's 



(11) 



Let us give a brief description of each of them 
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FIG. 4: Geometrical setting for {3Z, e) system 
t/j^^^: All cu's are chosen to be equal to ai, 

^(1) _ g-ai(ri+r2+r3) _ ^-|^2) 

It is a Heitler-London type function. This corresponds to coherent interaction between 
the electron and all centers. Supposedly, it describes the system at small intercenter 
distances and, probably, the equilibrium configuration. It is verified a posteriori. 

tp^"^^: Two a's are equal to zero and the remaining one is set to be equal to a2, 

It is a Hund-MuUiken type function. This function possibly describes the system at 
large distances, where essentially the electron interacts with only one center at a time 
thus realizing totally incoherent interaction. 

ip^^^: One a is equal to zero, two others are different from zero but equal to each other and 
to as, 

It is assumed that this function describes the system {2Z, e) plus center when a triangle 
is of a sufficiently small size. In fact, it is the Heitler-London function of (2Z, e)-system 
symmetrized over centers. 

ip^^^: One a is equal to zero and two others are different from each other being equal to q;4^5, 
respectively, 

^(4) _ g-a4ri-a6r2 _j_ g-Q;4T'2-a5n _j_ g-a4ri-a6r3 

_j_ g— a4r3— asn _j_ g— a4r2— as^'s _j_ g— o:4r3— a6r2 (15) 
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It is assumed that this function describes the system {2Z, e) plus one center. In fact, it 
is the Guillemin-Zener function of the {2Z, e)-system then symmetrized over centers. 
If 04 = as, the function ip^"^^ is reduced to ip^^\ If a4 = 0, the function ip^'^^ is reduced 
to tp^'^\ Hence, ^^^^ is a non-linear interpolation between -0*^^^ and 'ip^'^\ It has to 
describe intermediate intercenter distances. 

Two ct's are equal but the third one is different, 

^(5) _ ^-a(iri—a6r2-a7r3 _|_ ^—aeri-afir3-a7r2 
_j_ g-a6''2-«6»'i-a7''3 _j_ g-«6'''2-a6''3-a7'^l 

_j_ g— aera— «6ri-a7r2 _j_ ^—a(jrs-a6r2—a7ri (16) 

It describes a "mixed" state of three Z-hydrogen atoms. If ctg = ay, the function tp^^'' 
is reduced to tp^^\ If = 0, the function ip^^'^ is reduced to ip^'^\ If = 0, the 
function ip^^'^ is reduced to 'ip^^\ Hence, ip^^^ is a non-linear interpolation between ^p^^\ 
ip^'^^ and -0*^^). As function f|T5l) this is a type of Guillemin-Zener function and should 
describe intermediate intercenter distances. 

%p^^'>: All a's are different, 

^(6) _ ^—asri—agr2—aior:j _|_ ^—asri—agr:j—aior2 
^—asr2—agri—aior:j _|_ g— «8'''2— agr'a— aio'l 
-f- g~'^8^3~'^9^i~°i0''2 _j_ g-a8'"3-a9''2-aio''i (17) 

This is a general non-linear interpolation of all functions tp^^~^\ 

The total number of parameters of the function (fTOj) is equal to 15, where five are linear 
ones. Note that without a loss of generahty the parameter Aq in ffTOj) can be fixed, putting 
Aq = 1. We expect this function provides a relative accuracy ~ 10^3 in total energy. 

As a result of variational study for fixed Z the optimal geometric configuration is always 
the equilateral triangle. It was checked that this configuration is always stable with respect 
to small deviation. On Fig. [2] the total energy dependence for {3Z, e) at the equilibrium 
configuration is given. It is a smooth monotonous curve which ends at Z = Z^r < 1. At some 
charges this curve intersects with the energy curves for (2Z, e) and Z-atom, {Z,e). These 
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crossings separate domains of stability from different domains of metastability of the system 
(see below). On Fig. [3] the equilibrium distance between nearest static charges (the side of 
the equilateral triangle) is shown. It is a smooth curve which has a minimum i?™" = 2.413 
a.u. at Z = 0.391 and it grows to infinity with a decrease of Z. At critical charge Z = 0.9539 
the equilibrium distance, where the potential curve E = E{R) has the saddle point, is equal 
to 4.754 a.u. It is worth mentioning that the Req curves for {3Z, e) and {2Z, e) intersect at 
Z = 0.2670 with R^g = 2.506 a.u. 

It is interesting to study the approach of the total energy to the critical charge Z — > Zcr 
from below. In order to do it we use a general Puiseux expansion ([7]). Eventually, the 
behavior of the total energy close to critical charge Zcr, as a result of the fit, is given by the 
following terminated Puiseux expansion: 

E{Z) = -0.6954 + Q.2700{Zcr - Z) - 1.0357(Zc^ - Zf'^ - 1.3360(Zcr - Zf (18) 
-0.1350(Zcr - Zf^ + 2.3395(Zc^ - Zf - 1.8714(Zc^ - Zf'^ , 

where the critical point is 

Z'^^ = 0.9539 . (19) 

The fit (ITS]) is based on data from the domain Z G [0.80, 0.93] (19 points, see Table II). This 
behavior indicates that the critical point might be a square-root branch point. 

There are two points of crossing for the energy curve {3Z, e) at Fig. [2j The first one is 
Zcllss = 0.8269 with R^q = 3.234 a.u. for the crossing of {3Z, e) and the {Z, e). The second 

(2) 

one is ZcrLs = 0.5811 with R^q = 2.640 a.u. for the crossing of (3Z, e) and the {2Z,e) at 
Req = 2.008 a.u. For charges Z G (0.8269, 0.9537) for the triangular equilateral configuration 
the system is metastable with two decay channels 

(3Z,e) -^{Z,e) + Z + Z , 
(3Z,e) ^Z+ + Z , 

while for Z G (0.5811,0.8269) system is metastable with single decay channel 

(3Z,e) (2Z, e) + Z , 

and, finally, for Z < 0.5811 the system is stable. A study of the Puiseux expansions near 
Zcllss as well as Zcrlss from above and below show that they are the Taylor expansions which 
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z 




Fit 


0.10 


-0.057230 


-0.056880 


0.15 


-0.111367 


-0.111714 


0.20 


-0.173248 


-0.173581 


0.25 


-0.238656 


-0.238205 


0.30 


-0.304100 


-0.304235 


0.80 


-0.680137 


-0.680137 


0.82 


-0.682953 


-0.682953 


0.84 


-0.685194 


-0.685194 


0.86 


-0.686969 


-0.686969 


0.88 


-0.688417 


-0.688417 


0.90 


-0.689716 


-0.689716 


0.91 


-0.690385 


-0.690385 


0.92 


-0.691123 


-0.691123 


0.93 


-0.691991 


-0.691991 



TABLE II: Total energy Ej- of (3Z, e) in Ry at equilibrium vs Z obtained using (jlOp compared to 
the result of the fit (fT8|l . 

do coincide within the accuracy of data used and the obtained parameters of the fit. They 
do not give an indication that these points are branch points. It also is worth mentioning 
that the dependence of Rf,q on Z near ZlroL is also very smooth, see Fig. El 

Another question to rise is a behavior of the total energy near the critical point at Z = 
which is the singular point of the Schrodinger equation. Based on the fit of data from 
the domain Z G [0.1,0.7] (seven points, see Table II) we find that the Puiseux expansion 
becomes the Taylor expansion 

Et = -7.4257^2 + 19.3244Z^ - 19.4662Z^ + . . . . 

Such a behavior does not provide an indication to singular nature of the point Z = 0. 
However, the total energy can not be analytically continued to KeZ < 0. 
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C. (4Z,e) 



This Coulomb system consists of four static Z-charges and one electron, (4Z, e). It is 
worth anticipating that the most symmetrical configuration where Z— charges are placed on 
the vertexes of a tetrahedron, see Fig. is optimal. It was checked that this configuration is 
always stable with respect to small deviations. It is certain that for Z = 1 there is a bound 
state, the system Hl'^ does not exist. 




FIG. 5: Geometrical setting for (4Z, e) system 



Trial Functions. The variational method was used to obtain all numerical results. Trial 
function is taken in a form of linear superposition of three functions: 

3 



"Atrial = Ajtp 



(i) 



(20) 



where Aj are linear parameters. Each function ijA^'* is chosen in such a way as to describe 
different physical characteristics of the system |9|, [lO|. In general, ip^^^ has the form of a 
symmetrized product of four Is-Coulomb orbitals (Slater functions) 



_ g-oiri— a2r'2-03»'3-04r'4 



(21) 



Let us give a brief description of each of them: 
ip^^^ : All a's are chosen to equal to ai. 



(1) _ ^-ai{ri+r2+r-i+r4,) 
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(22) 



It is a Heitler-London type function. This corresponds to coherent interaction between 
the electron and all protons. Supposedly, it describes the system at small interproton 
distances and, probably, the equilibrium configuration. 

■0^^) : Three a's are equal to zero and the remaining one is set to be equal to a2, 

^(2) ^ g-aan ^ g-a2r2 ^ ^-a^rs ^ g-a2r4 _ (23) 

It is a Hund-MuUiken type function. This function possibly describes the system at 
large distances, where essentially the electron interacts with only one proton at a time 
thus realizing incoherent interaction. 

tp^^^ : All a's are different from each other, and different from zero. 
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(24) 



We can see that trial functions ip^^^ and ip'^'^^ are particular cases of the general trial function 

There might exist two critical charges which separates the domain of existence from non- 
existence of bound states. One such a critical charge is at Z = 0. Another one is at some 
Z — Zcr < 1. Calculations (see below) show that Zcr — 0.736 at Req — 6.50 a.u. , where 
the potential curve E — E{R) has the saddle point. Thus, the system (4Z, e) can exist 
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for charges < Z < Zcr- The energy dependence at equihbrium distance Req is a smooth 
function, see Fig. [2l The optimal geometric configuration is always the tetrahedron. It 
was checked that this configuration is always stable with respect to small deviations. It is 
a smooth monotonous curve which ends at Z = Z^ < 1. At some charges it intersects the 
energy curves for {3Z,e), {2Z,e) and Z-atom {Z,e). These crossings separate domains of 
stability from metastability of the system (see below). On Fig. [3] the equilibrium distance 
between nearest static charges (the side of the tetrahedron) is shown. It is a smooth curve 
which has a minimum Req = 2.485 a.u. at Z = 0.2218 and it grows to infinity with a decrease 
of Z. It is quite amusing that all three equilibrium distance curves for (4Z, e), (3Z, e) and 
(2Z, e) intersect for Z = 0.2670 with R^q = 2.506 a.u. 

Behavior of the energy as function of the charge close to critical charge Z^r is given by 
the terminated Puiseux expansion: 

E{Z) = -0.3368 - 0.2793(Z,^ - Z) - 1.5995(Ze^-Z)3/2 (25) 
+ 2.0214(Z„. - Zy + 0.9224(Zcr - Z)^/^ + . . . , 

where the critical point is 

= 0.736 . (26) 

The fit ( |25|) is based on data from the domain Z G [0.60,0.72] (12 points, see Table III). 
This behavior indicates that critical point might be a square-root branch point. 

There are three points of crossing for the energy curve (4Z, e) at Fig. [2l The first one is 
for the crossing of (4Z, e) and the {Z,e) energy curve at Zcllss = 0.6290 with Req = 4.187 
a.u. The second one is at zirlss = 0.4798 with Req = 3.08 a.u. for the crossing of (4Z, e) 
and the (2Z, e) energy curve at Req = 2.086 a.u. The third one is ZcrLs = 0.4065 with 
Req = 2.83 a.u. for the crossing of (4Z, e) and the {3Z, e) energy curve at Req = 2.413 a.u. 

For charges Z G (0.6290,0.7360) for the triangular equilateral configuration the system 
is metastable with three decay channels 

(4Z,e) {Z,e) + Z + Z + Z 
(4Z, e) ^ (2Z, e) + Z + Z 

(4Z, e) ^ (3Z, e) + Z (27) 
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z 




Fit 


0.10 


-0.092030 


-0.092063 


0.11 


-0.106530 


-0.106491 


0.12 


-0.121220 


-0.121190 


0.13 


-0.136020 


-0.136072 


0.14 


-0.151080 


-0.151070 


0.15 


-0.166150 


-0.166145 


0.60 


-0.411356 


-0.411387 


0.62 


-0.401055 


-0.401025 


0.64 


-0.389997 


-0.389987 


0.66 


-0.378426 


-0.378456 


0.68 


-0.366652 


-0.366675 


0.70 


-0.355022 


-0.354995 


0.72 


-0.344089 


-0.344021 



TABLE III: Total energy Et of (4Z, e) in Ry at equilibrium vs Z obtained variationally using the 
trial function ([20|) compared to the result of the fit ([25]) . 

For Z G (0.4798,0.6290) the system is metastable with two decay channels 

(4Z, e) ^ (2Z, e) + Z + Z , 

(4Z, e) ^ (3Z, e) + Z (28) 

For Z G (0.4065, 0.4798) system is metastable with one decay channel 

(4Z, e) ^ (3Z, e) + Z (29) 

and, finally, for < Z < 0.4065 the system gets stable. 

A question to rise is about behavior of the total energy near the critical point at Z = 
which is the singular point of the Schrodinger equation. Based on the fit of data from 
the domain Z G [0.,0.15] (seven points, see Table III) we find that the Puiseux expansion 
becomes the Taylor expansion 

Et = -14.3871^2 + 62.0529^3 - 102.4490Z^ + . . . . 
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Such a behavior does not provide an indication to singular nature of the point Z — 0. 
However, the total energy can not be analytically continued to ReZ < 0. 



IV. TWO-ELECTRON MOLECULAR SYSTEMS 



A. (2Z, e, e) 



The system {2Z, e, e) consists of two charged centers Z and two electrons. For Z = 1 it 



is the celebrated Ho molecule when for Z = 2 it is the Helium molecular ion He 



(++) 



which 



is metastable system. It is obvious that for large Z the system is unbound as well as for 
negative Z. Thus, there are two singular points: Z — where the potential "changes" sign 
and Zcr > 2 which is a critical point separating the domain of the existence from the domain 
of non-existence of the solutions in the Hilbert space. It seems natural that the ground state 
when exists is the spin-singlet state. Calculations (see below) show that the critical charge 
Zcr = 2.250 at Reg = 1.532 a.u. Thus, the system (2Z, e, e) exists for charge < Z < Zcr- 




FIG. 6: Geometrical setting for (2Z, e, e) system 



Trial Functions. To calculate the total energy of the {2Z, e, e) system as a function of 
charge Z the variational method is used. Exponential correlated trial functions with proper 
symmetrization are employed as well as their linear superposition. In general, the basic 
trial function ip^^^ is taken in the form of symmetrized product of four Is-Coulomb orbitals 
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(Slater functions) and correlation function in exponential form 12| 



+ e 



+ e 



+ e 



-Q:iria-oi2rib-a3r2a-air2b+iri2 



-a3ria-ct4rib-air2a-ct2r2b+iri2 



-a2ria-0!irib-a4r2a-0!3r2b+^ri2 



-a2ria-airib—a4r2a-a3r2b+iri2 



(30) 



Recently, it was shown that a linear superposition of three functions f l30p leads to the most 
accurate ground state energy for the H2-molecule among a few parametric trial functions. 




FIG. 7: Total energy Et in Ry vs Z for two-electron systems in equilibrium at i? = Req'. {2Z, e, e) 
(solid line) and {3Z,e,e) (dashed line), and for comparison for {Z,e,e) (long-dashed line). Both 
curves intersect at Z = 1.1767. Solid line ends at Z = ziV = 2.250. Dashed line ends at 



Z = Z, 



(2) 



1.433. 



With such a function the expected relative accuracy in total energy is ~ 10~^. The total 
energy E{Z, R = Req) is presented at Fig. [7] and the equilibrium distance is at Fig. |8l 
Both curves are smooth without any indication to a special feature at the physical charge 
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z 

FIG. 8: Equilibrium distance R^q in a.u. of system (2Z, e, e) (solid line) and (3Z, e, e) (dashed line) 
as functions of the charge Z. Curves cross twice at Z = 0.3460 and at Z = 0.6851. Solid line ends 
at Z = Z^r'^ = 2.250. Dashed line ends at Z = Z^? = 1.433. 

Z = 1,2. At charge Z = 1.1767 the energy curves for {2Z, e, e) and {3Z, e, e) intersect. The 
equihbrium distance is a smooth curve which has a minimum R^g = 1.264 a.u. at Z = 1.596 
and it grows to infinity with a decrease of Z. At critical charge Z = 2.250 the equihbrium 
distance, where the potential curve E = E{R) has the saddle point, is equal to 1.532 a.u. 

Behavior of the energy as function of the charge close to critical charge Z < Zcr is given 
by the terminated Puiseux expansion (see ([7])): 

E{Z) = -8.6835 + 5.5238(Zc^ - Z) - 0.2982(Zc^ - Zf/^ (31) 
- 0.3166{Z^r- Zf + 0.3577{Z^r- Zf/^ + ... , 

where the critical point is 

Z^r = 2.250 . (32) 
The fit ( 13T|) is based on data from the domain Z G [1.80,2.22] (12 points, see Table IV). 
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This behavior indicates that critical point might be a square-root branch point. 
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l.OU 




-O.oUOOl 


1 QD 
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-7.34768 


-7.34841 


2.10 


-7.87639 


-7.87623 


2.20 


-8.41058 


-8.41121 


2.22 


-8.52070 


-8.51953 



TABLE IV: Total energy Et of (2Z, e, e) in Ry at equilibrium distance vs Z obtained using the 
trial function (130]) compared to the result of the fit (f3T]) . 

There are two points of crossing for the energy curve (2^, e, e) displayed at Fig. [71 The 
first one is at ZcrLs = 1.7026 with R^q = 1.268 a.u. for the crossing of {2Z,e,e) and two 
atoms {Z,e). The second one is ZcrLs = 0.4501 with R^q = 2.126 a.u. for the crossing of 
{2Z,e,e) and the {2Z,e) at R^q = 2.119 a.u. and two atoms {Z,e). 

For Z G (1.7026,2.250) the system {2Z,e,e) is metastable, there is decay channel 

{2Z,e,e)^{Z,e) + {Z,e) , 

while for Z e (0.4501, 1.7026) system is stable and for Z < Zcross = 0.4501 it seemingly gets 
metastable again with two decay channels 

(2Z,e,e) ^ (2Z,e) + e , 

(2Z,e,e)-^(Z,e) + (Z,e) . 

About the last domain we are not certain due to insufficient accuracy of our calculations. 

A question to rise is a behavior of the total energy near the second critical point at Z = 
which is the singular point of the Schrodinger equation. Based on the fit of data from the 



22 



domain Z G [0.1, 0.6] (five points, see Table IV) we find that the Puiseux expansion becomes 
the Taylor expansion 



£^7 



-0.6533Z^ - 4.1162Z"^ + 2.5076Z^ + 



Such a behavior does not provide an indication to singular nature of the point Z = 0. 
However, the total energy can not be analytically continued to ReZ < 0. 

B. (3Z,e,e) 

The system (3Z, e, e) consists of three charged centers Z and two electrons. For Z = 1 
it is celebrated Hj]' molecular ion. It is obvious that for large Z the system is unbound as 
well as for negative Z. Thus, there are two singular points: Z = where the potential 
"changes" sign and Zcr > 2 which is a critical point separating the domain of the existence 
from the domain of non-existence of the solutions in the Hilbert space. It seems natural 
that the ground state when exists is the spin-singlet state. Calculations (see below) show 
that Zcr = 1.441 at Reg = 1.98 a.u. The optimal geometrical configuration at equilibrium is 
the equilateral triangle. It was checked that this configuration is always stable with respect 
to small deviations. Thus, the system {3Z, e, e) exists for charge < Z < Zcr- 




FIG. 9: Geometrical setting for (3Z, e, e) system 



Trial Functions. The variational method was used to obtain all numerical results. In 
general, the basic trial function tp^^^ has the form of a symmetrized product of six Is- 
Coulomb orbitals (Slater functions) and correlation function in exponential form (see 12 1 
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for a discussion), 
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(33) 



It is worth noting that a hnear superposition of three functions of a type ( 133|) leads to the 
most accurate energy for the Hj]'-molecule for lowest spin-triplet state in linear configuration 
^'Eu among a few parametric trial functions giving a relative accuracy ~ 10~^ 13|. 

With such a function f l5^ the expected relative accuracy in total energy is ~ 10~^. 
The total energy E{Z, R = Reg) of (3Z, e, e) is presented at Fig. [7] and the equilibrium 
distance is at Fig. [HI The optimal geometrical configuration at equilibrium is always the 
equilateral triangle. It was checked that this configuration is always stable with respect to 
small deviations. Both curves are smooth without any indication to a special feature at the 
physical charge Z = 1. At charge Z = 1.1767 the energy curves for {2Z,e,e) and (3Z, e,e) 
intersect. This crossing separates the domain of stability from metastability of the system 
(3Z, e, e): at Z > 1.1767 {3Z, e, e) can decay to (2Z, e, e) + Z. The equilibrium distance Reg 
is a smooth curve which has a minimum R^g = 1.643 a.u. at Z = 0.8981, it grows to infinity 
with a decrease of Z. Two equilibrium distances curves for (2Z, e, e) and {3Z, e, e) intersect 
twice for Z = 0.6851 with R^g = 1.6917 a.u. and for Z = 0.3460 with Reg = 2.4297 a.u. 

Behavior of the energy as function of the charge Z close to critical charge Z < Zer is 
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given by the Puiseux expansion (see ([7])): 

E{Z) = -3.6798 + 1.7613{Zcr - Z) - 0.5009(Zc^ - Zf^^ (34) 
+ 1.5164(Z,^ - + 0.6143(Zc^ - ^ . . . , 

where the critical point is 

Zcr = 1.441 . (35) 

The fit flMl) is based on data from the domain Z G [1.20, 1.42] (7 points, see Table V). This 
behavior indicates that critical point seems to be a square- root branch point. 
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-3.6268 
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-3.6436 


-3.6440 



TABLE V: Total energy Et of (3Z, e, e) in Ry at equilibrium vs Z obtained using ()33p compared 
to the result of the fit (pill . 

There are three points of crossing for the energy curve {3Z, e, e) displayed at Fig. [3 The 
first one is Zcllss = 1.3566 with Req = 1.881 a.u. for the crossing of {2Z,e) at Req = 2.406 
a.u. and the {Z, e). The second one is ZctIss = 1.3137 with R^q = 1.837 a.u. for the crossing 
with two atomic {Z,e) systems. The third one is ZcrLs = 1.1767 with Reg = 1.716 a.u. for 
the crossing with the (2Z, e, e) at R^q = 1.329 a.u. 
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For charges Z e (1.3566, 1.4407) the system is metastable with three decay channels 

(3Z, e, e) ^ (2Z, e) + (Z, e) 
(3Z, e, e) 2(Z, e) + Z 
(3Z, e, e) ^ (2Z, e, e) + Z 

For charges Z e (1.3137, 1.3566) the system is metastable with two decay channels 

(3Z, e, e) ^ 2(Z, e) + Z 
(3Z,e,e) {2Z,e,e) + Z 

For charges Z e (1.1767, 1.3137) the system is metastable with one decay channel 

(3Z,e,e) {2Z,e,e) + Z 

Eventually, for charges Z G (0.2989,1.1767) the system becomes stable. For charges Z < 
0.2989 the system can be either in the form (3Z, e) + e or (2Z, e) + {Z, e). The accuracy of 
our calculations do not allow us to make a definite statement. 

Another question to rise is a behavior of the total energy near the second critical point 
at Zcr — which is the singular point of the Schrodinger equation. Based on the fit of data 
from the domain Z e [0.1, 0.5] (five points, see Table V) we find that the Puiseux expansion 
becomes the Taylor expansion 

Et = -0.7198Z^- 11.9676Z^ + 14.0751Z^+... . 

Such a behavior does not provide an indication to singular nature of the point Z = 0. 
However, the total energy can not be analytically continued to ReZ < 0. 

Conclusions 

In this paper we calculated for the first time the critical charges of five simple 1-2 electron 
molecular systems: (2Z, e), (3Z, e), (4Z, e), (2Z, e, e), (3Z, e, e) under the assumption that 
the Z-charges are static and found their equilibrium configurations. It was also found 
that for all those systems the total energy and equilibrium distance vs Z are smooth curves 
without any indication to charge quantization. For all studied systems the optimal geometric 
configuration is the most symmetric being the equilateral triangle for (3Z, e), (3Z, e, e) and 
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tetrahedron for (4Z, e) . It was checked that this configuration is always stable with respect 
to small deviations. It seems natural to assume that for {AZ, e, e) the tetrahedron as the 
optimal geometrical configuration would occur. It would be interesting to study the optimal 
geometrical configuration for five (or more) Z-center cases, in particular, for (5Z, e). 

It is evident the existence of the critical charge for any one electron system {nZ, e), since 
the potential has a form V — —AZ + BZ^ with A,B>0 and always becomes negative if 
the charge Z is small enough. The critical charge behaves like Zcr,n oc with some a < 
at large n. However, the question about stability of {nZ, e) at Z < Zcr,n remains unclear to 
the present authors. Probably, a similar conclusion can be made for two electron systems. 
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